# Install if needed
#install.packages("tidyverse")
#install.packages("readr")
#install.packages("here")

library(tidyverse) # Versatile package for data analysis

library(readr) # read csv files
library(here) #easier to deal with relative paths 

library(ggthemes) # Themes for ggplot
library(directlabels) # Not sure, google it
library(ggrepel) # Add labels to plots

net_migration_map <- readRDS("datasets/net_migration_map.rds")
region.lab.data <- readRDS("datasets/region.lab.data.rds")

map_net_migration <- ggplot() + 

  geom_polygon(                                          # make a polygon
    data = net_migration_map,                                    # data frame
    aes(x = long, 
        y = lat, 
        group = group,          # coordinates, and group them by polygons
        fill = cat_mig),
    alpha=0.6
    ) 

map_net_migration

map_net_migration  + 

  geom_polygon(                                          # make a polygon
    data = net_migration_map,                                    # data frame
    aes(x = long, 
        y = lat, 
        group = group,          # coordinates, and group them by polygons
        fill = cat_mig),
    alpha=0.6
    ) + 

labs(title = "Net migration from 2013 to 2017",
     subtitle = " WHO Eastern Mediterranean Region",
     fill = "Net migration (thousands)", 
     caption = "Source: World Bank. Net migration is the net total of migrants during the period, less the annual number of emigrants, \n including both citizens and noncitizens. Data are five-year estimates. k- Thousands, m - Millions") +
  
  theme_map() 

map_net_migration  + 

  geom_polygon(                                          # make a polygon
    data = net_migration_map,                                    # data frame
    aes(x = long, 
        y = lat, 
        group = group,          
        fill = cat_mig),
    alpha=0.6
    ) + 

labs(title = "Net migration from 2013 to 2017",
     subtitle = " WHO Eastern Mediterranean Region",
     fill = "Net migration (thousands)", 
     caption = "Source: World Bank. Net migration is the net total of migrants during the period, less the annual number of emigrants, \n including both citizens and noncitizens. Data are five-year estimates. k- Thousands, m - Millions") +
  
  theme_map() +
  
  theme(plot.title=element_text( face="bold", 
                                 hjust=0.5,
                                 lineheight=1.2),
        plot.subtitle=element_text(hjust=0.5),
        legend.position = "bottom",
    panel.border = element_blank()
  )  +
  
  scale_fill_manual (
    values = c("#d7191c", "#fdae61", "#ffffbf", "#a6d96a", "#1a9641"),
    na.value = "gray55",
    name = "Net migration",
          guide = guide_legend(
            direction = "horizontal",
            keyheight = unit(2, units = "mm"),
            keywidth = unit(8, units = "mm"),
            title.position = 'top',
            # I shift the labels around, the should be placed 
            # exactly at the right end of each legend key
            title.hjust = 0.5,
            label.hjust = 1,
            nrow = 1,
            byrow = T,
            # also the guide needs to be reversed
            reverse = F,
            label.position = "bottom"
          )
  )

map_net_migration <- ggplot() + 

  geom_polygon(                                          # make a polygon
    data = net_migration_map,                                    # data frame
    aes(x = long, 
        y = lat, 
        group = group,          # coordinates, and group them by polygons
        fill = cat_mig),
    alpha=0.6
    ) + 

  
  geom_path(data = net_migration_map, 
            aes(x = long, 
                y = lat, 
                group = group), 
            color = "black", 
            size = 0.1) +

labs(title = "Net migration from 2013 to 2017",
     subtitle = " WHO Eastern Mediterranean Region",
     fill = "Net migration (thousands)", 
     caption = "Source: World Bank. Net migration is the net total of migrants during the period, less the annual number of emigrants, \n including both citizens and noncitizens. Data are five-year estimates. k- Thousands, m - Millions") +
  
  theme_map() +
  
  theme(plot.title=element_text( face="bold", 
                                 hjust=0.5,
                                 lineheight=1.2),
        plot.subtitle=element_text(hjust=0.5),
        legend.position = "bottom",
    panel.border = element_blank()
  )  +
  
  scale_fill_manual (
    values = c("#d7191c", "#fdae61", "#ffffbf", "#a6d96a", "#1a9641"),
    na.value = "gray55",
    name = "Net migration",
          guide = guide_legend(
            direction = "horizontal",
            keyheight = unit(2, units = "mm"),
            keywidth = unit(8, units = "mm"),
            title.position = 'top',
            # I shift the labels around, the should be placed 
            # exactly at the right end of each legend key
            title.hjust = 0.5,
            label.hjust = 1,
            nrow = 1,
            byrow = T,
            # also the guide needs to be reversed
            reverse = F,
            label.position = "bottom"
          )
  )

           
map_net_migration

map_net_migration <- ggplot() + 

  geom_polygon(                                          # make a polygon
    data = net_migration_map,                                    # data frame
    aes(x = long, 
        y = lat, 
        group = group,          # coordinates, and group them by polygons
        fill = cat_mig),
    alpha=0.6
    ) + 
  
  geom_text(data = region.lab.data,                                    # data frame
    aes(x = long, 
        y = lat, 
        label = region), 
            size = 2.5,
            hjust = 0.1,
    check_overlap = TRUE) +
  
  geom_path(data = net_migration_map, 
            aes(x = long, 
                y = lat, 
                group = group), 
            color = "black", 
            size = 0.1) +

labs(title = "Net migration from 2013 to 2017",
     subtitle = " WHO Eastern Mediterranean Region",
     fill = "Net migration (thousands)", 
     caption = "Source: World Bank. Net migration is the net total of migrants during the period, less the annual number of emigrants, \n including both citizens and noncitizens. Data are five-year estimates. k- Thousands, m - Millions") +
  
  theme_map() +
  
  theme(plot.title=element_text( face="bold", 
                                 hjust=0.5,
                                 lineheight=1.2),
        plot.subtitle=element_text(hjust=0.5),
        legend.position = "bottom",
    panel.border = element_blank()
  )  +
  
  scale_fill_manual (
    values = c("#d7191c", "#fdae61", "#ffffbf", "#a6d96a", "#1a9641"),
    na.value = "gray55",
    name = "Net migration",
          guide = guide_legend(
            direction = "horizontal",
            keyheight = unit(2, units = "mm"),
            keywidth = unit(8, units = "mm"),
            title.position = 'top',
            # I shift the labels around, the should be placed 
            # exactly at the right end of each legend key
            title.hjust = 0.5,
            label.hjust = 1,
            nrow = 1,
            byrow = T,
            # also the guide needs to be reversed
            reverse = F,
            label.position = "bottom"
          )
  )


ggsave( here::here("figures", 
                   "map_net_migration.png"), 
        dpi = "retina",
        width =15,            
        height = 10,          
        units = c("cm"),      
        )

library(osmdata)

available_tags("highway")
##  [1] "bridleway"              "bus guideway"           "bus stop"              
##  [4] "construction"           "corridor"               "crossing"              
##  [7] "cycleway"               "elevator"               "emergency access point"
## [10] "escape"                 "footway"                "give way"              
## [13] "living street"          "milestone"              "mini roundabout"       
## [16] "motorway"               "motorway junction"      "motorway link"         
## [19] "passing place"          "path"                   "pedestrian"            
## [22] "platform"               "primary"                "primary link"          
## [25] "proposed"               "raceway"                "residential"           
## [28] "rest area"              "road"                   "secondary"             
## [31] "secondary link"         "service"                "services"              
## [34] "speed camera"           "steps"                  "stop"                  
## [37] "street lamp"            "tertiary"               "tertiary link"         
## [40] "toll gantry"            "track"                  "traffic mirror"        
## [43] "traffic signals"        "trailhead"              "trunk"                 
## [46] "trunk link"             "turning circle"         "turning loop"          
## [49] "unclassified"
available_features()
##   [1] "4wd only"                "abandoned"              
##   [3] "abutters"                "access"                 
##   [5] "addr"                    "addr:city"              
##   [7] "addr:conscriptionnumber" "addr:country"           
##   [9] "addr:district"           "addr:flats"             
##  [11] "addr:full"               "addr:hamlet"            
##  [13] "addr:housename"          "addr:housenumber"       
##  [15] "addr:inclusion"          "addr:interpolation"     
##  [17] "addr:place"              "addr:postcode"          
##  [19] "addr:province"           "addr:state"             
##  [21] "addr:street"             "addr:subdistrict"       
##  [23] "addr:suburb"             "admin level"            
##  [25] "aerialway"               "aeroway"                
##  [27] "agricultural"            "alt name"               
##  [29] "amenity"                 "area"                   
##  [31] "atv"                     "backward"               
##  [33] "barrier"                 "basin"                  
##  [35] "bdouble"                 "bicycle"                
##  [37] "bicycle road"            "biergarten"             
##  [39] "boat"                    "border type"            
##  [41] "boundary"                "bridge"                 
##  [43] "building"                "building:fireproof"     
##  [45] "building:flats"          "building:levels"        
##  [47] "building:min level"      "building:soft storey"   
##  [49] "busway"                  "cables"                 
##  [51] "charge"                  "circuits"               
##  [53] "construction"            "covered"                
##  [55] "craft"                   "crossing"               
##  [57] "crossing:island"         "cuisine"                
##  [59] "cutting"                 "cycleway"               
##  [61] "denomination"            "diet"                   
##  [63] "direction"               "dispensing"             
##  [65] "disused"                 "disused:railway"        
##  [67] "disused:shop"            "drink"                  
##  [69] "drive in"                "drive through"          
##  [71] "driving side"            "ele"                    
##  [73] "electrified"             "embankment"             
##  [75] "embedded rails"          "emergency"              
##  [77] "end date"                "entrance"               
##  [79] "est width"               "fee"                    
##  [81] "fire object:type"        "fire operator"          
##  [83] "fire rank"               "foot"                   
##  [85] "ford"                    "forestry"               
##  [87] "forward"                 "frequency"              
##  [89] "fuel"                    "gauge"                  
##  [91] "generator:method"        "generator:output"       
##  [93] "generator:source"        "golf cart"              
##  [95] "goods"                   "hazmat"                 
##  [97] "healthcare"              "healthcare:speciality"  
##  [99] "height"                  "hgv"                    
## [101] "highway"                 "historic"               
## [103] "horse"                   "ice road"               
## [105] "incline"                 "industrial"             
## [107] "inline skates"           "inscription"            
## [109] "internet access"         "is in:city"             
## [111] "is in:country"           "junction"               
## [113] "landuse"                 "lanes"                  
## [115] "layer"                   "leaf cycle"             
## [117] "leaf type"               "leisure"                
## [119] "lhv"                     "line"                   
## [121] "lit"                     "location"               
## [123] "man made"                "maxheight"              
## [125] "maxlength"               "maxspeed"               
## [127] "maxstay"                 "maxweight"              
## [129] "maxwidth"                "military"               
## [131] "minspeed"                "mofa"                   
## [133] "moped"                   "motor vehicle"          
## [135] "motorboat"               "motorcar"               
## [137] "motorcycle"              "motorroad"              
## [139] "mountain pass"           "mtb scale"              
## [141] "mtb:description"         "mtb:scale:imba"         
## [143] "name"                    "narrow"                 
## [145] "natural"                 "noexit"                 
## [147] "non existent levels"     "note"                   
## [149] "nudism"                  "office"                 
## [151] "official name"           "old name"               
## [153] "oneway"                  "opening hours"          
## [155] "operator"                "organic"                
## [157] "oven"                    "overtaking"             
## [159] "parking:condition"       "parking:lane"           
## [161] "passing places"          "place"                  
## [163] "power"                   "produce"                
## [165] "proposed"                "protected area"         
## [167] "psv"                     "public transport"       
## [169] "railway"                 "railway:preserved"      
## [171] "railway:track ref"       "recycling type"         
## [173] "ref"                     "religion"               
## [175] "residential"             "resource"               
## [177] "roadtrain"               "route"                  
## [179] "sac scale"               "service"                
## [181] "service times"           "shelter type"           
## [183] "shop"                    "sidewalk"               
## [185] "site"                    "ski"                    
## [187] "smoothness"              "social facility"        
## [189] "start date"              "step count"             
## [191] "substation"              "surface"                
## [193] "tactile paving"          "tank"                   
## [195] "tidal"                   "toilets:wheelchair"     
## [197] "toll"                    "tourism"                
## [199] "tracks"                  "tracktype"              
## [201] "traffic calming"         "traffic sign"           
## [203] "trail visibility"        "tunnel"                 
## [205] "turn"                    "type"                   
## [207] "usage"                   "vehicle"                
## [209] "vending"                 "voltage"                
## [211] "water"                   "wheelchair"             
## [213] "wholesale"               "width"                  
## [215] "winter road"             "wires"                  
## [217] "wood"

getbb("Lisboa Portugal") # the coordinates of the city, save them 
##         min       max
## x -9.229836 -9.086333
## y 38.691399 38.796758


streets <- getbb("Lisboa Portugal")%>%
  opq()%>%
  add_osm_feature(key = "highway", 
                  value = c( "primary", 
                            "secondary", "tertiary")) %>%
  osmdata_sf()



small_streets <- getbb("Lisboa Portugal")%>%
  opq()%>%
  add_osm_feature(key = "highway", 
                  value = c("residential")) %>%
  osmdata_sf()

####
ggplot() +
  geom_sf(data = streets$osm_lines,
          inherit.aes = FALSE,
          color = "black",
          size = .4,
          alpha = .8) +
  geom_sf(data = small_streets$osm_lines,
          inherit.aes = FALSE,
          color = "black",
          size = .1,
          alpha = .2) +

  coord_sf(xlim = c(-9.296592, -8.976592), # Change the coordinates to your city 
           ylim = c(38.547751, 38.867751), # Change the coordinates to your city 
           expand = FALSE) +
  
  theme_void() 

ggsave( here::here("figures", 
                   "lisbon_map_white.png"), 
        dpi = "retina",
        width =15,            
        height = 10,          
        units = c("cm"),      
        )

indicators <- read_csv("datasets/indicators.csv")

indicators_sum <- indicators %>%
  filter( location_name %in% c("Afghanistan", 
                        "Bahrain", 
                        "Djibouti", 
                        "Egypt", 
                        "Iraq", 
                        "Iran",
                        "Jordan",
                        "Kuwait",
                        "Lebanon",
                        "Libya",
                        "Morocco",
                        "Oman",
                        "Pakistan",
                        "Palestine",
                        "Qatar",
                        "Saudi Arabia",
                        "Somalia",
                        "Sudan",
                        "Syria",
                        "Tunisia",
                        "United Arab Emirates",
                        "Yemen"
                        ),
         indicator_short=="Skilled Birth Attend",
         year_id<2019) 

indicators_sum <- indicators_sum %>%
  group_by(year_id) %>%
  mutate(WHO_EMR = mean(unscaled_value))
    
    
indicators_sum <- indicators_sum %>%    
    filter( location_name %in% c("Afghanistan", 
                        "Iraq", 
                        "Libya",
                        "Palestine",
                        "Somalia",
                        "Syria",
                        "Yemen",
                        "Sudan"
                        )
    ) %>%
  mutate(type = case_when( 
    location_name == "Afghanistan" ~ 0,
    location_name == "Iraq" ~ 1,
    location_name == "Libya" ~ 1,
    location_name == "Palestine" ~ 1,
    location_name == "Somalia" ~ 0,
    location_name == "Syria" ~ 1,
    location_name == "Yemen" ~ 0,
    location_name == "Sudan" ~ 0,
    ) )

plot.fit <- ggplot() +
  
      geom_line(data = indicators_sum,                 
               aes(x=year_id,               
                   y = unscaled_value,
                   group = location_name,
                   colour = factor(location_name)
#                   size=1 
               ), size=1 ) 
plot.fit 

plot.fit <- ggplot() +
  
      geom_line(data = indicators_sum,                 
               aes(x=year_id,               
                   y = unscaled_value,
                   group = location_name,
                   colour = factor(location_name)
#                   size=1 
               ), size=1 ) +
  
  scale_x_continuous(name="Year", limits=c(1990, 2020)) +
  
  labs(title = "Skilled birth attendance (%)",
     subtitle = "8 WHO Eastern Mediterranean Region states, 1990-2017",
     y = "Skilled birth attendance (%)",
     x = "Year", 
     fill = "States", 
     caption = "Source: IHME") 

plot.fit

plot.fit  +

  geom_vline(xintercept=1990,  color = "black") +
  
  theme_minimal() +
  theme(panel.grid.minor = element_blank(),
        panel.grid.major.y = element_blank(),
        axis.text = element_text(face = "bold", size = 8, color = "black"),
        axis.title = element_text(size = 10),
        plot.title = element_text(face = "bold", size = 12),
        legend.position = "none" # Play with it
        ) 

plot.fit  +
        geom_line(data = indicators_sum,                
               aes(x=year_id,                
                   y = WHO_EMR),
                   size=1,
               color="#ef6548",
               linetype = "dashed") + 
  
    geom_vline(xintercept=1990,  color = "black") +
  
  theme_minimal() +
  theme(panel.grid.minor = element_blank(),
        panel.grid.major.y = element_blank(),
        axis.text = element_text(face = "bold", size = 8, color = "black"),
        axis.title = element_text(size = 10),
        plot.title = element_text(face = "bold", size = 12),
        legend.position = "none" # Play with it
        ) 

plot.fit <- ggplot() +
  
      geom_line(data = indicators_sum,                 
               aes(x=year_id,               
                   y = unscaled_value,
                   group = location_name,
                   colour = factor(location_name)
#                   size=1 
               ), size=1 ) +
  
        geom_line(data = indicators_sum,                
               aes(x=year_id,                
                   y = WHO_EMR),
                   size=1,
               color="#ef6548",
               linetype = "dashed") +
  
      geom_text_repel(data  = subset(indicators_sum, 
                                    year_id == 2018),
                               aes(x=year_id,               
                   y = unscaled_value,
                       label = location_name),
                    nudge_x = 2,
                    size=3,
#                  box.padding   = 0.35, 
#                  point.padding = 0.5,
                  segment.color = 'grey50') +
  
  scale_x_continuous(name="Year", limits=c(1990, 2020)) +
  
  labs(title = "Skilled birth attendance (%)",
     subtitle = "8 WHO Eastern Mediterranean Region states, 1990-2017",
     y = "Skilled birth attendance (%)",
     x = "Year", 
     fill = "States", 
     caption = "Source: IHME") +
  
  annotate(
    geom = "curve", 
    x = 2002, 
    y =     58, 
    xend = 2004, 
    yend =  60, 
    curvature = .3, 
    arrow = arrow(length = unit(2, "mm"))) +
  annotate(geom = "text", 
           x = 2000, 
           y = 60, 
           label = "Iraq conflict", 
           hjust = "left",
           size=3) + 
  annotate(
    geom = "curve", 
    x = 1998, 
    y =     82, 
    xend = 1995, 
    yend =  82, 
    curvature = .3, 
    arrow = arrow(length = unit(2, "mm"))) +
  
    annotate(geom = "text", 
           x = 1996, 
           y = 80, 
           label = "Palestine conflict", 
           hjust = "left",
           size=3)  +
  
      annotate(geom = "text", 
           x = 2014, 
           y = 86.5, 
           label = "WHO-EMR", 
           hjust = "left",
           size=2.8)  +
  
  
  geom_vline(xintercept=1990,  color = "black") +
  
  theme_minimal() +
  theme(panel.grid.minor = element_blank(),
        panel.grid.major.y = element_blank(),
        axis.text = element_text(face = "bold", size = 8, color = "black"),
        axis.title = element_text(size = 10),
        plot.title = element_text(face = "bold", size = 12),
        legend.position = "none"
        ) 

plot.fit 

ggsave( here::here("figures", 
                   "skills_rate.png"), 
        dpi = "retina",
        width =16,            
        height = 10,          
        units = c("cm"),      
        )

who_conflict <- read_csv("datasets/who_conflict.csv")

who_conflict <- who_conflict %>%
  filter(metric_name=="Number")

who_conflict <- who_conflict %>%
  group_by(year) %>%
  mutate(sum_deaths= sum(val),
         percent=(val/sum_deaths)*100
         )

library(viridis)
cbp1 <- c("#eff3ff", "#fb6a4a", "#bdd7e7", "#6baed6",
          "#3182bd", "grey")


who_conflict %>% 
  ggplot( aes(x=year, 
              y=percent, 
              fill=location_name)
          ) +
    geom_area( )  

library(viridis)
cbp1 <- c("#eff3ff", "#fb6a4a", "#bdd7e7", "#6baed6",
          "#3182bd", "grey")


who_conflict %>% 
  ggplot( aes(x=year, 
              y=percent, 
              fill=location_name)
          ) +
  
    geom_area()  + 
  
  scale_fill_manual(values=cbp1) +
  
  labs(title = "Percent of total deaths by year from war and terrorism from 1990 - 2017",
       subtitle = "World Health Organization regions",
       y = "Percent of total deaths from war and terrorism",
       x = "Year", 
       fill = "", 
       caption = "Source: IHME") 

library(viridis)
cbp1 <- c("#eff3ff", "#fb6a4a", "#bdd7e7", "#6baed6",
          "#3182bd", "grey")


who_conflict %>% 
  ggplot( aes(x=year, 
              y=percent, 
              fill=location_name)
          ) +
  
    geom_area()  + 
  
  scale_fill_manual(values=cbp1) +
  
  labs(title = "Percent of total deaths by year from war and terrorism from 1990 - 2017",
       subtitle = "World Health Organization regions",
       y = "Percent of total deaths from war and terrorism",
       x = "Year", 
       fill = "", 
       caption = "Source: IHME") +
  
    annotate("text", 
           x = 2012, 
           y = 40, 
           label = "Eastern Mediterranean Region ", 
           color = "black", 
           size=3) +
  
    annotate("text", 
           x = 1998, 
           y = 75, 
           label = "African Region", 
           color = "black", 
           size=3) +
  
      annotate("text", 
           x = 2005, 
           y = 10, 
           label = "South-East Asia Region", 
           color = "black", 
           size=3) +

  theme_minimal() +
  
  theme(panel.grid.minor = element_blank(),
        panel.grid.major.y = element_blank(),
        axis.text = element_text(face = "bold", size = 8, color = "black"),
        axis.title = element_text(size = 10),
        plot.title = element_text(face = "bold", size = 12),
        legend.position = "bottom"
        ) 

ggsave( here::here("figures", 
                   "area_chart.png"), 
        dpi = "retina",
        width =16,            
        height = 12,          
        units = c("cm"),      
        )